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ABSTRACT 

We use a suite of high-resolution N-body simulations to study the properties, abun- 
dance and clustering of high mass halos at high redshift, including their mass assembly 
histories and mergers. We find that the analytic form which best fits the abundance of 
halos depends sensitively on the assumed definition of halo mass, with common defi- 
nitions of halo mass differing by a factor of two for these low concentration, massive 
halos. A significant number of massive halos are undergoing rapid mass accretion, with 
major merger activity being common. We compare the mergers and mass accretion 
histories to the extended Press-Schechter formalism. 

We consider how major merger induced star formation or black hole accretion 
may change the distribution of photon production from collapsed halos, and hence 
reionization, using some simplified examples. In all of these, the photon distribution 
for a halo of a given mass acquires a large scatter. If rare, high mass halos contribute 
significantly to the photon production rates, the scatter in photon production rate can 
translate into additional scatter in the sizes of ionized bubbles. 



1 INTRODUCTION 



Observations of the anisotropy of the cosmic microwave 
background (CMB) radiation have given us unprecedented 
knowledge of the very early Universe and dramatically con- 
firmed the picture of large-scale structure as arising from 
the gravitational amplification of small perturbations in a 
Universe with a si gnificant cold dark matter component 
l|Smoot et al.l [l993 ). In this model the ionization history 
of the Universe has two main events, a 'recombination' at 
z ~ 10^ in which it went from ionized to neutral and a 'reion- 
ization' during z ~ 7 — 12 in which the radiation from early 
generations of collapsed objects was able to ionize the inter- 
galactic medium. The former event is strongly constrained 
by the CMB. A new generation of instruments will soon 
allow us to probe this second event: "the end of the dark 



estimates of the number and properties of virialized dark 
matter halos at high redshift, the focus of this paper. At 
z — 10, halos with M > 10^ Mq are expected to be bi- 
ased similarly to very massive clusters (M > lO'^^ h~^MQ) 
today, with the most massive and recently formed halos 
growing rapidly and merging frequently. We explore some 
properties of these collapsed halos at a high redshift using 
a suite of high resolution, coUisionless, N-body simulations. 
We pay particular attention to merger rates and mass ac- 
cretion histories with an eye to applications for reionization. 
We also co mpare the N-body results w ith the predictions of 
the oft-used IPress fc Schechterl l| 19741 ) formalism. 

If halo mergers are accompanied by a temporary in- 
crease in photon production (due eithe r to starbursts 
increased black hole accre t ion e.g. [ Carlberd 199C ; 



20011: ICoorav & BartonI 


20061: iFan, CariUi & Keatinsll2006l: 


Furlanetto. Oh & BriEES 


20061). 



Barnes fc Hernquis tI [l99ll. 1 19961: iMihos fc Hernguisd 



■ 1994 

W96; Kauffmann fc Haehneltl 12000: Cav aliere fc Vittorinil 
20001) we expect reionization to be influenced by the merger 



Since at reionization a very small fraction of the mass 
affected each and every baryon in the Universe, reionization 
is particularly sensitive to the distribution and behavior of 
collapsed structure. We expect that the ionizing sources are 
situated in large (Tvir > 10*K or M > W h'^MQ) dark 
matter halos where the gas can cool efficiently to form stari|3- 
Models for the sources of reionization thus often start with 



^ We will only consider Pop II stars here; Pop III stars, which can 
form in the absence of m etals in smaller halos, arc expected to be 
less likely by redshift 10 llYoshida. Bromm fc Hernguist 2004) . 



and accretion history of dark matter halos, beyond just the 
fact that more massive halos emit more photons. With a sim- 
ple model of star formation we show that merger- induced 
scatter in photon production may be significant, with the 
production rates acquiring a substantial tail to large pho- 
ton production rates. Since the massive halos are relatively 
rare, this individual halo scatter is expected to translate into 
a scatter of photon production rates inside ionized regions, 
changing the bubble distribution. 

The outline of the paper is as follows. In ^we describe 
the N-body simulations. The basic halo properties are de- 
scribed in ^ along with the results for mergers and mass 
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Figure 1. Illustrative merger trees for two halos with masses of 1.2 (left) and 2.9 X 10^° H-'^Mq (right). Time runs upwards in steps of 
10 Myr, from z = 12.7 (bottom) to z = 10 (top) and the age of the Universe (in Myr) is shown at every second step. At each time the 
area of the symbol is proportional to the halo mass, with masses decreasing to the right in each group, and lines show the progenitor 
relationship. The leftmost branch shows the main trunk of the tree. The halo at left has a (major) 1:2 merger at the last time step, while 
the main trunk of halo at right has a 1:2 merger at the first time step, a 1:6 two steps later and then only smaller mergers after that. 



gains and the comparison to Press-Schechter. The conse- 
quences of this merging in a simple model for photon pro- 
duction are elucidated in ^and we summarize and conclude 
in 31 



2 SIMULATIONS AND PARAMETERS 

We base our conclusions on 5 dark matter only N-body sim- 
ulations of a ACDM cosmology with Qm = 0.25, JIa — 0.75, 
h — 0.72, n — 0.97 and as = 0.8, in agreement with a 
wide array of observations. The initial conditions were gener- 
ated at z = 300 using the Zel'dovich approximation applied 
to a regular, Cartesian grid of particles. Our two highest 
resolution simulations employed 800'^ equal mass particles 
{M = 2 X 10® and 1.7 x lO^/i-^M©) in boxes of side 25 
and 50/i~^Mpc with Plummer equivalent smoothings of 1.1 
and 2.2/i.^^kpc. They were evolved to z = 10 using the 
TreePM code described in IWhitd (120021) (f or a comparison 
with other codes see iHeitmann et al.ll2007l ). We ran 3 ad- 
ditional, smaller simulations in a 20/i~^Mpc box, one with 
600^ particles and two with 300'^ particles (each started at 
z = 200). A comparison of the boxes allows us to check for 
finite volume, finite mass and finite force resolution effects. 
We shall comment on each where appropriate. 

The phase space data for the particles were dumped 
at 15 outputs spaced by 10 Myr from z — 12.7 to 2: = 10 
for all but the largest box. The lower resolution of the 
largest box makes it less useful for merger trees, so it was 
sampled for only subset of these output times, ending at 
z — 10. For each output we generate a cata log of halos 
using the Friends-of-Friends (FoF) algorithm (iDavis et al.l 
[l985) with a linking length, b, of 0.168 times the mean inter- 
particle spacing. This partitions the particles into equiv- 
alence classes, by linking together all particle pairs sepa- 
rated by less than b. The halos correspond roughly to par- 
ticles with p > 3/(27r6^) ~ 100 times the background den- 
sity. We also made catalogs using a linking length of 0.2 
times the mean inter-particle spacing, which we shall dis- 
cuss further below. We found that the FoF algorithm with a 
larger linking length had a tendency to link together halos 



whi ch we would, by eye, have characterize d as separate (see 
also iDavis et al.lll985l : ICole fc Lacevll 19961 . for similar discus- 
sion). This problem is mitigated with our more conservative 
choice of b. 

For each halo we compute a number of properties, in- 
cluding the potential well depth, peak circular velocity, the 
position of the most bound particle (which we take to define 
the halo center) and Afigo, the mass interior to a radius, 
riso, within which the mean densi ty is 18 times the back- 
ground densitfl As discussed in Fwhitej (^001, 2002) and 
iHu fc Kravtsovri|2003l ). the choice of halo mass is problem- 
atic and ultimately one of convention. We shall return to 
this issue in the next section. 

Merger trees are computed from the set of halo catalogs 
by identifying for each halo a "child" at a later time. The 
child is defined as that halo which contains, at the later time 
step, more than half of the particles in the parent halo at 
the earlier time step (weighting each particle equally). For 
the purposes of tracking halos this simple linkage between 
outputs suffices (note that we do not attempt to track sub- 
halos within larger halos, which generally requires greater 
sophistication). Two examples of the halo merger trees are 
given in Fig. [TJ where we see a rich set of behaviours, in- 
cluding major and minor mergers and many body mergers. 
From the merger trees it is straightforward to compute the 
time when a halo 'falls in' to a larger halo, the number and 
masses of the progenitors etc. 

Due to finite computational resources, all N-body sim- 
ulations must trade-off computational volume for mass res- 
olution. By running multiple simulations we can overcome 
this to some extent, but not entirely. We have chosen to 
slightly under-resolve the low mass (Tvir ~ WK) halos m 
order to simulate a slightly larger volume, since our focus 

^ Note this is simply a definition of halo mass, not the halo finder. 
We still use FoF particles to define the group centers. However 
given the center we use all of the particles in the simulation when 
determining Migo- Our Migo masses should thus be comparable 
to the sum of the particles in an 50(180) group - a common def- 
inition that employs both the 50(180) halo finder and definition 
of mass. 
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2.5 3 3.5 4 4.5 no" 10^ IQi" 



iv= 6^ (t)/a(M) M,3„ or M,„, (M^/h) 



Figure 2. The peak height, u - which governs the abundance, 
clustering and merging behavior in analytic models - for z = 0, 
4, 7 and 10. For example, objects with u = 3 have Af ~ 4 X 
10* h-^MQ at 2 = 10 but Af ~ 6 X 10^'' H-^Mq at z = 0. 



will be on the more massive halos which have more frequent 
major mergers. Under reasonable assumptions (see below) 
between i — | off all photon production occurs in halos 
more massive than 10^ h~^MQ at z = 10, and we easily re- 
solve these objects with the 25 h~^Mpc simulation which we 
use for the bulk of the paper. 



3 HALO PROPERTIES 

3.1 Halo abundance and clustering 

The highest mass objects in our volume have mass ~ 
10^°/i~^Mq and radii of several tens of kpc. At z = 10 
these halos are analogous to rich clusters today, being re- 
cently formed and rare: Fig. [2] shows the mass as a function 
of peak height, i' = (5c/cr(M), at 2 = 10, 7, 4 and 0. The 
threshold Sc{t) is defined as 1.686/D{t), where D is the lin- 
ear growth factor normalized to unity at z = and a^{M) 
is the variance of the mass computed using linear theory at 
z = 0. In our cosmology 5c{z = 10) ~ 13.8. Due to the flat- 
ness of the dimensionless power on the scales of interest, the 
slightly red initial spectrum and the low clustering ampli- 
tude, the characteristic mass, M*, where a ~ 5c, is 0{1)Mq 
at z = 10, so all of the halos we consider are 3> M*. 

One of the most basic and useful quantities we 
can derive from the simulations is the mass function, 
the spatial abundance of halos as a function of mass. 
High redshift mass functions have been stud ie d by many 
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■■2oo7r Heitmann et al' 

Trac & Ccn 2006; Ilicv ot al. 2005, 2006a; Maio et al 



Zahn et all I2OO7I ; iLukic et al.l 120071 ) and iLukic et al 



1 2OO7I ) offer a comprehensive summary of recent work. Most 
previou s work finds mas s func ti ons which are better fit 
by the [S heth fc Tormeiil (|l999l '). Ijenkins et al] l|200ll ') or 



^^^etal t2Q0j) = W^d that the appropr^e mass 
function to use depends primarily on the definition of mass 
chosen and definitions which at z ~ give very similar mass 
functions can give quite different ones sX z = 10. 



Fi gure 3. The m ass functions fo r our box compared to the 
War ren et all ||2005|) (dashed l ii ie). Ijenkins et al.l l200lh (dash- 
dotted line), ISheth fc TormenI l ll999l ) (dotted hne) and Press- 
Schechter (solid lower line) mass functions for our cosmology. 
N-body Afiso results are plotted as solid symbols for the 800^ 
(25/i-lMpc, squares, SO/i-^Mpc, hexagons) and 600^ (20/1-^ 
Mpc, triangles) runs. Only masses where there are more than 
10 halos in the box and where resolution effects are unimportant 
are shown. Open symbols denote the analogous FoF(0.2) mass 
functions for the same simulations. 



We show the mass function(s) from our three high- 
est resolution simulations in Fig. [S] If we use as our mass 
estimator the sum of the particle masses in the FoF(0.2) 
gro ups (open symbol s ) the n w e find good agreeme nt with 
the ISheth fc TormenI (|l99i) or I Jenkins et al.l (|200ll ) forms. 
This is the procedure followed by most of the groups abov^l. 
However if we choose instead to use Migo as our mass es- 
timator (filled symbols) we find a different mass function. 
Although this mass functio n shows a mark ed excess o f high 
mass halos compared to the |Pres~fc Schechterl (| 19741 ) form, 
it is a better fit than the alternate forms mentioned above. 
Agreeably, for the scales plotted, the M180 mass function is 
independent of the initial FoF group catalog used to define 
the centers about which Migo is determined. This is not too 
surprising as the group centers hardly change and the num- 
ber of "small" groups which split off of larger FoF groups 
as the linking length is decreased is tiny compared to the 
number of low-mass "field" halos. The differences in mass 
functions then comes primarily from the definition of the 
masses of the found objects. Comparing halo by halo the 
FoF(0.2) masses are almost twice M180, though the differ- 
en ce depends o n ma ss. A similar difference was also noted 
bv lReed et all (|2007l ) as a shift to lower abundance at fixed 
mass when comparing an FoF(0.2)-based mass function to 
that of a different halo finder. We believe the primary is- 



^ Due to the finite size of our boxes the mass function is slightly 
suppressed at high mass. We can estimate this suppression using 
(extended) Press-Schechter theory, assuming we have simulated 
the conditional mass function within a region of exactly mean den- 
sity on the mass scale of the box. The mass functions plotted have 
been corrected for this expected suppression, which ranges from 
< 1% t o 22% over the m ass range plotted. See, e.g. iLukic et al.l 
ll2007l)MReed et al.l ||2007^ for further discussion. 
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Figure 4. Density profiles of the 5 most massive halos in the 
800^ run a.t z = 10. The masses range from 1 — 3 X lO^'' h ^Mq. 
The 2 halos with the flatter profiles (short and long dashed lines) 
correspond to the S"^"* and 4*'' most massive halos and both have 
had a major merger (greater than 1;6) within the previous 10 Myr. 
Halo sizes (rigo) are below 100 /i~^kpc for all of the halos shown. 
The solid line, offset, shows an isothermal sphere profile (p oc r~^) 
for comparison. This indicates why FoF(0.2) masses assuming 
an isothermal profile may be expected to disagree with SO(180) 
masses as discussed in the text. 
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Figure 5. (Top) The halo-dark- matter cross correlation, ihmix) 
(top), for halos with (comoving) number density lO"'^", 10~^'^, 
10~^" and lO""'^ /i^Mpc"^ (open symbols from top to bot- 
tom) from our 50h~^Mpc simulation. The solid squares show 
the dark matter correlation function, ^mm(''). The ratio, b(r) = 
£.hm{^)/imm{r), is shown in the lower panel. 



sue is not the halo finder, but the mass definition. Their 
second halo finder assigns masses which are essentially our 
Mi8o. The mass discrepancy is much larger for these halos 
at z; = 10 than it is for group and cluster- sized halos at the 
present day (e.g. Figure 11 in IWhitdl2002l ). 

The mass differences are quite interesting. The his- 
torical argument for choosing FoF(0.2) was that the FoF 
group finder selects particles approximately within a den- 
sity 3/(27r&^) ~ 60 times the mean density. For a singu- 
lar isothermal sphere profile (p cx r~^) and a critical den- 
sity Universe the mean enclosed density is thus ISOpcrit, 
in accord with argument s based on spherical top-hat col- 
lapse (e.g. iPeacoc M ll998l ). At z = 10 the Universe is close 
to critical density so we might expect the FoF(0.2) and 
Miso mass functions to agree better than at lower z where 
180p ~ 45pcrit. However, we are focusing on very high mass 
halos which have only recently formed at z = 10. They are 
therefore less centrally concentratecfl than a 'typical' halo. 
This can also be seen in Fig. 2] where halo profiles are com- 
pared to the isothermal sphere profile. As the halo profiles 
are less steep than the isothermal sphere form assumed in 
the argument above, this leads to the differences in mass 
between the FoF(0.2) and 5'O(180) definitions. 

By contrast, we find that the FoF(0.168) mass function 
is very similar to the Migo points plotted, and a halo by 
halo comparison shows that the two masses agree to within 
20-30 per cent. As we go down the mass function, to more 
concentrated, low er mass halos, we ex pect FoF(0.2) to better 
match Afi8o (e.g. lCole fc Lacevl[l996l ). 

In general, given the strong dependence of the mass 



^ They corr espond roughly to c :^ 2 — 5 f or halos of the form 
proposed bv lNavarro. Frenk fc White! lll997t) . 



function upon the mass definition, and the ambiguity in 
this quantity in many analytic treatments, significant care 
must be taken when making predictions for the abundance 
of halos. Even if we decide to treat all halos as a simple 
1-parameter family, it is likely preferable to make compar- 
isons with some quantity more directly related to observ- 
ables (such as circular velocity, halo virial temperature or 
potential well depth) or to discuss statistics as a function of 
number density rather than mass. 

Like rich clusters we expect that these massive ha- 
los, in the process of formation, will not lie on the usual 
'vacuum' virial relation 2KE=PE, where KE and PE re- 
fer to the potential and kinetic energy respectively. In 
fact we find that 2KE/PE~ 1.4 for halos in the range 
10* — 10^° h~^ M(^, very similar to the value found for rich 
clusters to day (iKnebe fc MulleillT999l : ICohn fc WhTt^bOOSl : 
ih aw et al. Il2007ir A similar 'excess' k inetic energy was also 
found by Ijang-Condell fc Hernguis^ HooJ) for lower mass 
halos. The ratio is larger than unity because of the stead y 
accretion of material onto the cluster l|Cole fc Lacevlll996l ). 

Fig. \S\ shows the clustering of the dark matter and the 
halos from our 50/i"^Mpc run. We plot the auto-correlation 
function of the dark matter and the cross-correlation of the 
halo centers with the dark matter respectively. The latter is 
both less subject to nois«[3 from our small sample of massive 
halos and more applicable to understanding how radiation 
from the halos would infiuence the surrounding mass. The 



^ There is essentially no shot-noise for the dark matter, ^mm(f), 
and jackknife errors on the cross-correlation, ^hmy a^re a few per- 
cent for the samples shown. Jackknife drastically underestimates 
the errors from finite volume however. 
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Figure 6. A comparison of the large-scale bias measured for the 
mass thresholded samples of Fig. \5\ with a number of theoret- 
ical models: the bias of the P rcss-Schcchtcr mass function (as 
computed by lEfstathiou et "ah! 1988; Cole & Kaiser 1989, solid), 
the Sheth-Tormen mass fu nction dSheth fc Tormen| il999i. dashed) 
and the fitting function of ISheth. Mo fc TormenI (|200ll . dotted). 
Although the mass fun ction is in good agreement with that of 
ISheth fc TormenI l|l99^ ). their bias formula underestimates the 
clustering of the rarest halos. 



ratio of the cross- to auto-correlation functions defines tiie 
scale dependent bias, bh{r). 

The mass auto-correlation function is in good agree- 
ment between the 25 and 50 ft^^Mpc boxes up to 1 /i~^Mpc, 
with ^ from the 25ft.~^Mpc box falling below that of the 
50/i~^Mpc box beyond this scale. The 20/i~^Mpc box has 
noticeably less power over a wide range of scales. For the 
masses where we can compare and for the range of linear 
scales plotted, the halo-mass cross-correlation functions of 
the 25 and 50/i~^Mpc boxes are in excellent agreement, so 
we have shown the results only for the 50 h~^Mpc box. 

Our halo samples are mass thresholded, however by us- 
ing number density as our marker we largely avoid the is- 
sues of mass definition discussed earlier. The differences in 
bias at fixed n for the different mass choices, arising from 
the scatter between different mass definitions, is only a few 
percent. Taking 6^(1.5 /i~^Mpc) as the asymptotic value, 
the la rge-scale bias is in good a greement with the mod - 
els of iPress fc Schechte3 Il974l): lEfatathiou etlll (1 1988 



ICole fc Kaiseil ||1989| ): IMo fc White! (Il996l):ljiiid lll998h and 
~ 30% higher than tha t of 'Shcth fc Tormen' ^199^). Those 
of ISheth. Mo fc Tormen (2001) and Tinker ct al. (2001) lie 



in between. (The model of ISeliak fc WarrenI |200j) only ex- 
tends up to 100 times the non-linear mass, where 6 ~ 3, and 
it not applicable to our results.) To make contact with the 
earlier literature we plot in Fig. [6] the bias as a function o f 
peak h eight, z/, obtained from n using the Sheth fc Tormen 
mass function. When computing b(> v) in the sim- 
ulation we rank order the halos by FoF(0.2) mass in order 
to best match the chosen mass function. This mass function 
is then used when analytically computing the halo-weighted 
bias 6(> u) from each o f the analytic forms wh ich provide 
h{v). Because of this the lSheth fc TormenI l| 19991 ) bias func- 
tion is the only one which would give an average bias of 
unity when integrated over u. 
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Figure 7. (Top) The fraction of halos with M > 10^ H'^Mq 
which have had a 1:10, 1:5 or 1:3 merger (top to bottom) back to 
the lookback time shown in all 4 of our simulations. (Bottom) The 
fraction of halos with M > 10^ h~^MQ which have a large mass 
gain (mj /rrii > 1.1, 1.2, 1.33) vs. time. Here we show only the 
two highest resolution simulations for clarity. Both plots would 
coincide if all mergers were 2-body within 10 Myr. 



Similar trends for rare halos to have larger bias than 
the mode rn fits predict have been seen at lower red- 
shift (e.g. IShen et all l2008l: IWhite. Martini fc CohnI l2008l : 



lAngulo. Baugh fc LacevI 20081 . for recent work) but we must 
also remember that bh{l-5 h~^Mpc) is likely higher than 
bh{r — > oo) so the degree of overshoot is hard to quantify 
precisely. As expected, the clustering strength is an increas- 
ing function of mas s (|Kaiseij Il984l : lEfstathiou et al.l Il988l : 
ICole fc Kaiseill 19891) . or a decreasing function of halo abun- 
dance. 



3.2 Mergers and Mass Gains 

We now consider the hierarchical assembly of the dark mat- 
ter halos through merging and accretion. We shall use the 
800'', 25 h~^Mpc simulation since it provides both high mass 
resolution and a representative volume. Since our progeni- 
tor relationships are based on particles in the FoF groups, 
we use the FoF(0.168) masses for consistency. As discussed 
earlier, for our massive halos these masses are within 20 — 30 
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per cent of Miso and none of our conclusions depend sensi- 
tively on this choice. Fig. [7] shows the fraction of halos with 
10^ < M < 10^° h~^MQ which have experienced at least 
one major merger as function of lookback time, in intervals 
of lOMyrs. We show three different definitions of 'major' 
merger, where the largest two progenitors of the halo have 
ratios below 1:10, 1:5 or 1:3. Mergers are frequent but not 
ubiquitous - not all halos have had a major merger within 
140 Myrs, but many have. The fraction decreases for smaller 
mass ratios and for lower mass halos, as expected. 

We can also consider mass gains between time steps, 
often denoted in the literature as mf/rrii where rrii is the 
mass of the largest progenitor at the earlier time and ruf 
is the mass of the halo under consideration. Mass gains are 
sometimes used as a proxy for mergers. Fig. [7] shows those 
halos whose mass increased by at least 10, 20 or 33 per cent 
as a function of lookback time. The top and bottom panels 
of Fig. [7] would be identical if all mergers were two body and 
there was no smooth accretion. As can be seen in Fig.[T]this 
is not the case; Fig. [7] quantifies this difference for major 
mergers. 

The Press-Schechter model predicts the evolution of the 
mass function, and it can be extended to make predictions 
for the time history of halos. This "excursi on set formalism" 
is often called extended Press-Schechter (iBond et al.lll99l|: 
Bower! Il99ll : iLacev fc Coll Il993l. Il994l ; iKitavama fc Sutol 
1996h and denoted EPS - see IZentneiri|2006h for a recent 
review. Although it is analytically tractable, it has many 
inconsistencies and does no t compare particularly well to 
N-bo dy simulations (see e. g.lSheth fc Pitman! 1997 : TormerJ 
1998 : Somorvillo ct al. 20g^j7 



Cohn. Bada fc Whitel l20ol 



Benson. Kamionkowski fc Hassani 20051 : Li et al.ll2006l ). For 
example, Li et al. ( 20061 ) found that with EPS halos of mass 



10'^ - W^h'^MQ at z ~ formed later thaii in N -body 
simulations (but see IPercival. Miller fc Peacock|[2000l . for a 
slightly different quantity). In Fig. [8] we compare the N- 
body mass accretion histories for massive halos at z = 10 to 
a model by I Miller et al.l l|2006l ) based on EPS which predicts 
almost exponential growt h with redshift. (Other analytic 
models also exist, see e.g. iNeistein. van den Bosch fc Dekel 
l|2006t ) for a summary and comparison, there are some dis- 
crepancies between these which are not yet fully under- 
stood.) We find that EPS predicts mass growth which is 
too rapid also for the high mass, high redshift regime stud- 
ied here. The N-body mass accretion histories are relatively 
well fit b y an exponenti a l in z - a growth model also pro- 
posed bv lWechsler et al.l l|2002l ) on the basis of N-body sim- 
ulations of g alaxy-sized halos at low z - but the coefficient 
predicted bv lMiller et"al] l|2006l ) is larger than measured in 
the simulations. 

Perhaps the most common use of EPS is to predict 
merger rates, and EPS has been used in this context in 
several recent models of reionization. To compare the EPS 
predictions with our simulations we computed merger rates 
using only our last (lOMyr) time step, taking for any halo 
with z = IQ mass within M[ to M/ -I- AAf/ the distribution 
of progenitor^, Mprog- The EPS prediction can be found in 
the Appendix. 
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Figure 8. The mass accretion history for halos in the range 
(5 - 8) X 10** /i-^Mq from the SOO^ (d ashed) and 600^ (d otted) 
simulations and the functional form of iMiller et al. I 1 I2OO6I . solid) 
based on EPS. 
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Figure 9. EPS (solid) and simulation (dotted and dashed) results 
for the number of progenitors of halos with Mf = (4 — 4.5) X 
10* h~^MQ as a function of Mprog- 



We show a representative example of N{A'[prog)/N{A4f) 
for Mf in the range (4-4.5) x 10® H'^Mq in Fig.|9] For most 
of the range the agreement is reasonably good. At the low 
mass end EPS significantly underpredicts the number of pre- 
decessors found in our simulations (see also iPerc ival 2001j). 
At the high mass end the EPS rate starts to climb rapidly, 
eventually diverging unphysically. These trends are indepen- 
dent of the final mass chosen, or the definition of mass used. 
The EPS formula as progenitor mass goes to zero also di- 
verges, which we could not approach due to our finite mass 
resolution, but the mass weighted EPS calculation is finite 
at both endfl There is another notable difference between 
EPS and our simulation. Though it is relatively small, our 
time step is still too large for all mergers to be truly 2-body 
(see Fig. [T|, as implicitly assumed by EPS. A large fraction 



° We thank D. Holz for suggesting this as a useful comparison 
quantity. 



^ We thank Jun Zhang for emphasizing this. 
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(20-50 per cent, depending on Mf) of the halos are actually 
produced in 3 (or more)-body mergers. 

Finally we also looked for evidence that recently merged 
halos clustered differently than randomly chosen halos of 
the same mass. The correlation function of 1:2 or 1:3 
mergers appeared to be slightly (< 10%) enhanced at 
1 Mpc compared to the random sample, but the num- 
ber of merged halos was too small for this to be statisti- 
cally meaningful. The effect thus appears to be modest, if 
present at all, j ust as was found for lower redshift, high- 
mass halos (e.g. iGottlober et all I2OO2I: iPercival et all [20031 : 



IScannapieco fc Thackerl 20031 : IWetzel et all bOOTh . This sug- 
gests that the clustering of massive halos does not depend 
strongly upon their recent merger history. This in turn sig- 
nificantly eases the modeling of merger-related processes, 
such as enhanced photon production during reionization 
which we now discuss. 



4 REIONIZATION EFFECTS 

The rate of photon production in a galaxy can be 
enhanced by mergers, which can trigger starbursts or 



possibly ac cretion 



onto a 



black hole which 



ma^ 



be 



present (e.g. Carlberd 1990: Barnes & Hernquist 1991, 1996"; 



Mihos fc Hernguistl Tl994. .1996: .Kauffmann fc Haehnelt, 



200d : ICavaliere fc Vittorinill2000D . It is reasonable to antici- 



pate that the mergers of large dark matter halos could have 
similar effects on the photon production rate of the sources 
within them. We will make this assumption, and then con- 
sider the consequences of the merger rates computed above 
for the photon production distribution. 

We frame our discussion in terms of a sim- 
ple but promising model for reio n ization proposed by 
iFurlanetto. Zaldarriaga fc HernquistI (|2004lJ ). though our 
result is true more generally. In these models, a halo of a 
given mass m (in units of some reference mass) is consid- 
ered a source of photons with rate 



— = (t(m)m 



(1) 



Usually (t is taken to be mass independent, scale 



,2/3 



or transition from m 



2/3 



to m" at M 



10^" /i"^Mq (Furlanctto, McQ uinn fc Hernquisj|2005l . mo- 
tivated bv iKauffmann et al.ll2003 ). A region around these 
halos is taken to be ionized if the photons within it are 
sufficient to ionize all the interior mass. Some extensions 
also give recombinations spatial and/or temporal depen- 
dence and inc orporate this into findi ng the bubble properties 
ijFurlanetto fc Oh 2005; Fu rlanetto. McQuinn fc HernquistI 
l2005l : ICohn fc Chanei,2007.'). or incorporate Eg. (ffl) into N- 



body simulation s (|niev et al.1 l2006al : iMcQuinn et all l2007l : 
IZahn et al.ll200'7h . Under these assumptions the morphology 
of ionized regions can be computed from the photon produc- 
tion rate and spatial distribution of dark matter halos. 

A first step at including halo mergers within the 
above formalism (and its generalizations) was presented in 
ICohn fc Chan^(|2007^ . Those calculations were based on the 
Press-Schechter formalism, and so could only provide aver- 
age numbers of mergers for halos in a given mass range; 
scatter was computed by assuming that the mergers had a 
Poisson distribution. With our simulations we are able to 



check these assumptions and significantly extend this work 
because we have access to the detailed merger history of 
each halo. This allows us to go beyond their analytic esti- 
mates to explicitly calculate the full distribution of photon 
production for a halo of mass m, taking into account the 
distribution of histories and their associated (and difi^erent) 
photon production rates for a fixed m. 

From the merger tree for each halo at 2; = 10 (t = tobs) 
we identify which progenitors had at least one major merger 
(greater than 1:3 or 1:10), and the time fmorgc they occurred. 
We include all of the mergers in the tree and we place fmorgo 
at random within the lOMyr interval between the relevant 
outputs. Each of these mergers is allowed to contribute "ex- 
cess" photons beyond those which would automatically be 
assigned to the halo on the basis of its z = 10 mass, Alh, but 
the number of photons contributed is exponentially attenu- 
ated with an e-folding time r. The "excess" photon produc- 
tion is thus proportional to 



E 

merge 



s)/r 



(2) 



where the sum is over all halos which have undergone a 
major merger and we take a = 1 or 5/3. The exponential 
decay is motivated by modeling of starbursts, e.g. IConselicd 
hence the subscript s. We also consider another vari- 
ant, including all halos with major mergers within t of iobs, 
with no attenuation: 



Mt 



E 



Me{tobs-tn 



(3) 



where 0{x) = 1 if a; > 0, 1/2 if a:: = and zero otherwise. 
We denote this by a subscript bh, to indicate photon pro- 
duction by black holes, which might have their photon pro- 
duction rate increase over time and then decay once the fuel 
is exhausted. Assuming a step-like function is a crude first 
approximation to this uncertain physics. In all cases we take 
the quiescent photon production to depend on the z = 10 
halo mass with the same index, a, as Ms. We note this pre- 
scription might cause some overcounting if many mergers 
occur within a short time period and the gas becomes de- 
pleted from the earliest ones. A more refined model would 
account for the evolving baryon budget within the halo, but 
our treatment is sufficient for the purpose of illustration. 

The relative amplitudes of these two modes of photon 
production dep e nd on a number of different factors (see e.g. 
ICohn fc ChaQ (|2007l l for discussion and summary of esti- 
mates at these redshifts) but a factor /3 ~ 5 is not unreason- 
able for starbursts and could be even larger for black holes. 
The total photon production is thus enhanced by a factor 



(4) 



for the "starburst" prescription, or its analogue Ms Mbh 
for the "black hole" prescription. In principle both can con- 
tribute. We considered the two effects separately, their com- 
bination is straightforward. 

Figure [10] shows a typical example of the cumulative 
distribution of enhancement factors, Eq. Q. We took the 
starburst form, 1:3 mergers, a = 5/3, r = 75 Myr and l3 = 5, 
but other cases are very similar. The enhancement distribu- 
tion is extended, with a long tail to high £mrg and a peak at 
those halos which have not merged. About half of the halos 
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Figure 10. The cumulative probability for enhancement of pho- 
ton production over the quiescent case, £mrg, for one toy model. 
In this example, mergers (1:3) enhance the quiescent photon pro- 
duction rate (~ m^^^) at ^ = 10 with /3 = 5 and r = 75Myr 
(see text). Roughly 80 per cent of these halos have some en- 
hancement, shown for 2 bins in halo mass (with the range in 
logiQ{M/h~^ Mq) shown in each panel). There are 1195 and 136 
halos in the low- and high-mass bins, respectively. The peak at 
Emrg = 1 comes from the 20 per cent of the halos which have 
had no 1:3 mergers in the previous 140 Myr. Trends for other 
parameters and parameterizations are discussed in the text. 



tons. Even choosing oc m", such halos contribute ~10 
per cent of the photons. The number density of such halos 
is ~ 0.03 /i^Mpc"^. Bubble radii in different models range 
over several orders of magnitude. A middle-of-the-road esti- 
mate is 3/i~^Mpc, which would contain about 3 halos with 
M > 10^ Mq. The bubble radius would also be larger 
than the correlation length of our halos, so clustering is only 
expected to change this number by a factor of order unity. 
A small number of halos contributing a large fraction of 
the photons means that scatter in their photon production 
should affect the properties of the bubbles. 

Our calculation is relatively crude, but it suggests that 
the inclusion of mergers into a more refined model of reion- 
ization could alter the distribution of ioni zed regions. For 
models based on approximate dynamics (e.g. McQuinn et al.l 
l20Q7l : lMesinger fc Furlanettdl2007l : [Zahn et al.ll2007l ). a pos- 
sible first step would be to assign a merger history to the 
sources at random. This is accurate to the extent that re- 
cently merged halos are not spatially biased with respect 
to a random sample of halos of the same mass. For mod- 
els which marry the analytic model to dark matter sim- 
ulatiouilfl the merger history is known, so only the pho- 
ton production rate needs to be modified. More complex 
simulations involving radiative transfer will need to follow 
the photon production history as the halo s evolve, per- 
haps using a semi-analytic model ( s uch as in Benson, et al.l 
I2OOII : ICiardi. Stoehr fc Whitell2003l : iBenson. et al.ll2006D . A 
full-blown simulation including radiative transfer and N- 
bo dy in a large enough volume is still out of reach (but 
seelSo kasian et al.ll2003l: iKohler. Gnedin fc HamihonI I2OO5I : 
iTracfc Cen..2006l ; llliev et aLlbOOGal lbT for recent progress^ . 



have twice the photon production, while 20 per cent have no 
enhancement. Choosing a larger (3 increases the size of the 
enhancement, but does not qualitatively change the form 
of the distribution. Similarly, changing a or r changes the 
detailed form of the distribution but not its character. Ha- 
los down to 10* h~^MQ show a very similar distribution of 
enhancements. By contrast the model for black hole accre- 
tion produces a bimodal distribution, as the "early" mergers 
contribute relatively more than in the case of the starbursts, 
leading to a second peak. 

Even though the scatter in photon contributions from 
halo to halo is large for a given mass, if a large number 
of such halos are found in a bubble, their contributions to 
the photon numbers will tend to the mean, allowing the 
distribution to be replaced by the average. Precisely count- 
ing the number of halos of a given mass and the com- 
bined photon scatter inside a typical bubble is unfortunately 
self-referential: changing the ionization properties (including 
mergers) changes the bubble sizes and thus the number of 
halos within. Different assumptions about the nature of the 
sources and their feedback can give drastically different bub- 
ble sizes, a nd the relative import a nce of high vs. low er mass 
halos (e.g. iMcQuinn et al.1 [2OO7I : IZahn et al.1 l2007h . Given 
these uncertainties we consider properties in an average vol- 
ume, for illustration. 

For quiescent photon production and C, cx rn^^^ , ana- 
lytic estimates such as Press-Schechter give that halos with 
M > 10^ h~^MQ contribute between i — | of all pho- 



5 SUMMARY AND CONCLUSIONS 

Using 5 N-body simulations with different sized boxes and 
particle loads we considered the abundance, clustering and 
assembly histories of high mass halos at high redshift. We 
present results specifically for z = 10, but the evolution of 
the populations is smooth and the results will be similar 
at slightly higher and lower redshift. Like the halos of rich 
groups or clusters today the halos we consider are in the 
process of forming, growing rapidly through accretions and 
mergers. We found that they had larger velocity dispersions 
than a naive application of the virial theorem would pre- 
dict, due to a surface pressure from infalling material. Being 
recently formed, the halos were not very centrally concen- 
trated, leading to a factor of two difference between FoF(0.2) 
masses and Afiso. When measured against Miso we found 
our halo abundances were closer to the Press-Sche chter fit- 
ting formula than that of ISheth fc TormenI l| 19991 ). though 
the simulations had more high mass halos than the analytic 
form. If FoF(0.2) masses we re used instea d, the mass func- 
tion approached that of .Sheth fc TormenI (|l999.), in agree- 
ment with earlier work. This discrepancy indicates that an- 
alytic models which assign an observable to halos of a given 
size need to pay particular attention to the marker of halo 
size employed. 



° Unfortunately our simulation volumes are too small to provide 
converged answers for this step with the existing runs. 
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The high mass halos were significantly clustered, and we 
calculated the halo bias by taking the halo-mass cross cor- 
relation and dividing by the matter auto-correlation func- 
tion. Our rare halos were more clu s tered than the recent 
models of Sheth fc TormenI (Il999l ): ISheth. Mo fc TormenI 
ll200lll: [Tinker et all (|2005|) and closer to t he models of 
lEfsta thiou et alJ (l l988l ): ICole fc Kaise3 l| 19891 ) : [Mo fc White! 
||1996); Jing (1998). 

Merging is common, though not ubiquitous, in high 
mass halos ai z — 10. Major mergers, with progenitor mass 
ratios less than 1:3, occurred within 140 Myr of 2 = 10 
for more than half of halos with M > W^h'^MQ. We 
looked at the fraction of halos undergoing mergers for a va- 
riety of lookback times and progenitor ratios, finding more 
mergers for more massive halos, longer lookback times or 
less extreme merger events. Mass gains, parameterized by 
mf/rrii, showed similar trends even though not all merger 
events were two body within our 10 Myr time step. The EPS 
model provides a reasonable description of the progenitor 
mass distribution, though it underpredicts the number of 
low mass progenitors and diverges as Mprog Mh- The 
mg iss accretions histories predicted by EPS, as calculated 
by iMiller etHI |2006i l. provide only a qualitative guide to 
the mean mass accretion histories seen in our simulations. 

At z = 10 reionization is expected to be underway 
due to photon production from astrophysical objects which 
formed in collapsed halos. Within the context of a simple 
model which associates mergers with an increase in photon 
production rate the photon production distribution devel- 
oped a high rate tail due to recently merged halos. Includ- 
ing these photon production enhancements will likely drive 
scatter in photon production at fixed halo mass. Since the 
number of massive halos within a typical ionized bubble can 
be small, this scatter in photon production could well trans- 
late into additional scatter in bubble sizes and it would be 
very interesting to include this efi^ect in approximate models 
of reionization. If the recently merged halos are not spatially 
biased with respect to other halos of the same mass, includ- 
ing these effects in models, even those without merger trees, 
should be straightforward. 

JDC would like to thank D. Holz and S. Koushiappas 
for conversations. We thank O. Zahn for helpful comments 
on an early draft, and S. Furlanetto, CP. Ma, L. Miller, 
E. Neistein and J. Zhang for comments on the final draft. 
The simulations were analyzed on the supercomputers at 
the National Energy Research Scientific Computing center. 
MW was supported by NASA. We thank the referee for en- 
couraging us to run the 50/i~^Mpc box. 



APPENDIX 

We compared the progenitor mass distribution from the 
simulations with the predictions of the extended Press- 
Schechter formalism. To compute the distribution in this 
formalism we need the number of progenitors with mn < 
rui < mi2 at time t — At for a given range of final masses, 
m/i < nif < m/2, at time t. The EPS theory gives 



where the number of hal os of mass M, N(M,t), i s given by 
the usual mass function (|Press fc Schechteij[l974 ): 



N {mprog, At) 
At 



dM N{M,t) 



X {M/Mi)Pi{Mi~* M,t) 



dMi 



N{M, t) 



1 po 5c{t) 
2^M o^{M) 



da-' 



dM 



exp 



2a2(M) 



(6) 



with Po the background density. To get Eq. ((5)1 the number 
of halos at a given time is multiplied by the fraction of halos 
that have mass M ai t and have jumped from mass Mi 
within 5t: 



{M/Mi)Pi{Mi M-t)dMi3t = 



(7) 



(5) 



where Act^ = a{M{)'^ - a{M)^. To produce the EPS predic- 
tions discussed in ^3.2\ we take St — At = 10 Myr, assuming 
that this time is small. 
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